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We present a new approach to study the thermodynamic properties of d-dimensional classical systems by 
reducing the problem to the computation of ground state properties of a d-dimensional quantum model. This 
classical-to-quantum mapping allows us to deal with standard optimization methods, such as simulated and 
quantum annealing, on an equal basis. Consequently, we extend the quantum annealing method to simulate 
classical systems at finite temperatures. Using the adiabatic theorem of quantum mechanics, we derive the rates 
to assure convergence to the optimal thermodynamic state. For simulated and quantum annealing, we obtain the 
asymptotic rates of T(t) ~ (pN)/(kB log£) and ~f(t) ~ (Nt)~ c ^ N , for the temperature and magnetic field, 
respectively. Other annealing strategies, as well as their potential speed-up, are also discussed. 
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An outstanding issue in combinatorial optimization is the 
classification of problems according to their computational 
complexity. Typically, one defines a cost function that needs 
to be minimized and the question is how the number of re- 
sources (e.g., time) to determine the minimum scales with 
the problem size N. Long time ago it has been recognized 
that certain physics problems can be cast in this language. 
For example, it has been shown that the computation of the 
ground state energy (or the partition function) of classical 
three-dimensional spin glasses belongs to the class of NP- 
complete problems UJ], i.e. there is no known algorithm that 
can find the solution with polynomial (in N) resources. After 
all, the number of possible microscopic configurations of the 
system increases exponentially with the system size N and, 
unless certain symmetries reduce the complexity, one has to 
search in an exponentially large state space. This simplifica- 
tion happens, for example, in the two-dimensional Ising spin 
glass 1 2] (or any planar graph or lattice). 

Simulated annealing (SA) j^l and quantum annealing 
(QA) |4[ represent general algorithmic strategies to attack 
these optimization problems. The basic idea consists in find- 
ing the solution to the optimization problem as a limit of an 
effective physical process which uses additional variables or 
dimensions, and where the cost function is identified with a 
Hamiltonian H of a classical physical system. In SA one 
introduces temperature T as a tunable parameter: Initially 
the system is heated and next T(t) is slowly decreased to- 
wards zero, eventually converging to the ground (lowest en- 
ergy) state, whose energy equals the cost function. In QA, 
however, a time-dependent ad-hoc external magnetic field of 
magnitude j(t) is added to H, such that the total Hamiltonian 
can be interpreted as that of a quantum system. The (quantum) 
annealing process consists of slowly decreasing j(t) from a 
large value towards zero, while keeping T = 0. Since a 
quantum system in d dimensions can be mapped onto another 
classical system in d + 1 fSj, effectively in QA one is adding 
one extra space-dimension to the problem. In both strategies, 



the annealing procedure is essential to converge to the desired 
(ground) state, as it avoids getting stuck in local minima, using 
less resources than other optimization methods | 6]. 

In this Letter, we propose new algorithms to study the ther- 
modynamic properties of classical systems (including frus- 
trated systems, such as spin glasses). The crux of the method 
consists of mapping the classical d-dimensional problem into 
a quantum problem of the same dimensionality, and then us- 
ing techniques similar to those of QA to solve the latter. Our 
particular mapping allows us to unify the methods of SA and 
QA, and extend them to: i) study arbitrary classical models 
at T > and ii) study new annealing schemes. From this 
classical-quantum mapping perspective, any annealing strat- 
egy differs by the choice of path in (quantum) Hamiltonian 
space. Computation of thermodynamic properties of the clas- 
sical model amounts then to computation of ground state prop- 
erties of the mapped quantum model. Our approach can be 
readily implemented on a classical computer (CC) by using 
existent stochastic methods, such as Green's Function Monte 
Carlo 01, or by simulating the corresponding time-dependent 
Schrodinger equation. Since the proposed algorithms are 
based on a slow change of interactions in the quantum sys- 
tem, the rate at which these can be changed to assure conver- 
gence to the desired final state is determined by the adiabatic 
theorem |8]. Remarkably, we will show that for the path cor- 
responding to SA, the adiabatic condition yields to the result 
obtained by Geman and Geman on the rate of convergence to 
the optimal (ground) state of the classical system 0. 

For simplicity, we study classical models defined on a lat- 
tice (or graph), where a variable <jj = ±1 is defined on each 
site (vertex) j, and is related with the states of a physical spin- 
1/2. Any spin configuration of the 2 N possible ones is denoted 
as [a] = [<7i, • • • , crjv], where N is the total number of sites 
(or problem size). An energy functional E[a] (cost function) 
is defined on the lattice and its value depends on the state [a]. 
For example, in the Ising model, E[<r] = J^ij JijO~iO-j, where 
two interacting spins i and j contribute (— Jy) to the en- 



2 



ergy if they are in the same (different) state(s). In the canoni- 
cal ensemble, the expectation value of a thermodynamic vari- 
able A at temperature T is given by 



Z(T) 



(1) 



M 



where Z(T) = e is the partition function and (3 = 

M 

1/(/cbT), with ks the Boltzmann's constant. 

Any classical (finite-dimensional) spin model on a lattice 
can be associated with a quantum one, defined on the same 
lattice, by mapping every classical state [a] into a quantum 
state | [a]}. In this way, the energy functional maps into a 
Hamiltonian operator H. For spin- 1/2 models, H is given 
by mapping <jj — ► u\ in E[<r], where a\ is the Pauli opera- 
tor acting on the jth site. For example, H = Jij^Wi m 
the Ising model. The TV-spin (unnormalized) quantum state 
\ip(T)) = e- fm l 2 J2 [a] IM) (i.e., the Gibbs state), satisfies 



(A) = tr [pA] = 



(V(T)|i^(T)) _ 



(2) 



where p 



-0H 



/ Z(T). The operator A is determined by 



mapping the thermodynamic variable A, as described above. 
Then, [A, H] = 0. 

The state \tp(T)) can be shown to be the ground state of a 
family of quantum Hamiltonians H q (T) 111(111 . which are de- 
fined on the same lattice. Each of these Hamiltonians can be 
connected through a similarity transformation to a possible 
transition matrix M q (T) of a Markovian process leading to 
the thermal distribution: H q {T) = 1- e-> 3H / 2 M q (T)e l3H / 2 . 
Interestingly, the interactions appearing in H q (T) are of com- 
parable range to the interactions of the classical model. There- 
fore, a finite T phase transition of a d— dimensional classical 
system can then be identified with a quantum phase transi- 
tion of a d— dimensional quantum model. Thus, constructing 
a specific H q (T) and studying its ground state properties is 
of paramount importance as it will allow us to build different, 
yet more efficient algorithms to determine the thermodynamic 
properties of the classical system. 

We obtain the simplest form of H q (T) in the following 
way. First, notice that the Pauli operator a 3 x (i.e., the spin- 
flip operator acting on the jth site) satisfies erje~' 9ii ' 2 <7;j. = 
e PH je -0H/2 ( yj g [i ) jv]. The Hamiltonian Hj contains the 
terms in H having the operator a\ (i.e., the terms in H that an- 
ticommute with a j x ). Moreover, aiJ2[a] N) = E[ CT ] IM)> 
and H((T)\tp(T)) = 0, where H((T) = a{ - e^. In 
the basis determined by the states |[er]), the off-diagonal el- 
ements of H 3 q (T) are non-negative, and the coefficients ap- 
pearin g in \%j){T)) are all positive. The Perron-Frobenius the- 
orem lllll guarantees then that for T > 0, \ip{T)) is the 
unique ground state of the irreducible quantum Hamiltonian 
H q (T) = -xT,j HJ q( T )- The coefficient x = e~ 0p , with 
p ?s maxj \Hj\ = 0(1), is set for normalization purposes in 
order to satisfy \H q (T — > 0)| < oo. At this point, we would 



like to emphasize the simplicity of our particular mapping: 
The thermodynamic properties of any spin- 1/2 classical sys- 
tem can be obtained by studying the ground state properties of 
a spin- 1/2 quantum model, with classical interactions deter- 
mined by T and H (i.e., the classical system), and an external 
(homogeneous) transverse field of magnitude \- Remarkably, 
this field generates quantum fluctuations that are in one-to- 
one correspondence with the classical fluctuations at temper- 
ature T. In particular, H q (T -> oo) w (N - J2j a l)' so its 
ground state has all spins aligned along the external field, i.e. 
\ip(T — > oo)) ?s \ i a })- This quantum state can be iden- 
tified with the completely mixed state in the classical model. 



In the limit of low T we obtain H q (T ~ 0) f=a ^E, e 



whose expectation value is minimized by the ground state(s) 
of the classical model, i.e. \4>(T ~ 0)) is also a lowest energy 
state of H. 

To illustrate these results, we consider the homogeneous 
one-dimensional Ising model H = J Eilj <T z <7 z +1 - m 
this case, W q (T) = o{ - x 2 - xy{<j{- y ai + o{al +1 ) - 
y 2 al" 1 ai +1 , with x = cosh(/3J), and y = sinh(/?J). The 
Hamiltonian H q (T) denotes then a frustrated quantum Ising 
model, with next-nearest-neighbor interactions, and a trans- 
verse magnetic field of magnitude \- Such a frustration for- 
bids the existence of an ordered quantum phase unless T — * 0. 
This result is related to the non-existence of an ordered phase 
at finite temperature in the classical model. 

Within our context, we can interpret the SA procedure as a 
(real time) quantum evolution where we start from the initial 
quantum state ^{T — > oo)) « S[<r] IM)' an d next we c ' e " 
crease the interaction parameter T(t) (related to the tempera- 
ture of the classical model) in H q (T). If such an evolution is 
performed adiabatically, we remain in the desired ground state 
\i/j(T(t))) at any time t. Therefore, the gap A(T) between the 
ground and first excited states of H q (T) plays an important 
role on the rate at which T(t) must be decreased. This gap 
can be shown to satisfy A(T) > 2V2^Ne- ( - fJ P+ 1 ^ N = A(T). 
Such a lower bound can be determined using the inequalities 
in Ref. O and considering that (N - H q (T)) N is a strictly 
positive operator. It is based on the worst-case scenario (i.e., 
for the most general form of H), so it is expected to be im- 
proved depending on the nature of the interactions of the clas- 
sical system, such as translational invariance. The rate of the 
evolution is then determined by the adiabatic condition 1 8 ] 



W m (T(t))\&rH q (T)\il>(T(t))) 



A 2 n (T(t))^WW) 



- d t T = e, < t < T, 
(3) 

where e determines an upper bound to the probability of 
finding the system in any other (normalized) excited eigen- 
state \ij] m (T)) of H q (T), A m (T) is the energy gap between 
\ip m (T)) and \^{T)) (e.g., A X (T) = A(T)), and T is the to- 
tal time of the evolution. The Ihs of Eq. (0 can be bounded 
above by pN[2k B T 2 A(T)]- 1 \d t T\. To see this, note that 

d T H q {T)\i>{T)) = [d T (-f3H/2),H q (T))\i;(T)), (4) 
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as —f3H/2 generates the translations of \ip(T)). Therefore, 
\(^ m (T)\d T H q (T)\^(T))\ _ \(i> m (T)\H\ip(T))\ 



mechanics. To show this, we define a quantum Hamiltonian 



A m (T) 



(2fc s T2) 



(5) 



with \(ip m (T)\H\ip(T)}\ < pNyJZ{T). This upper bound is 
not necessarily tight. Equation (0 implies a resource require- 
ment of T ps 0[l/eA(T)] instead of T ps C[l/eA 2 (T)], 
which is the common resource scaling associated with an adi- 
abatic evolution. [Nevertheless, both scalings will yield to 
similar asymptotic behavior for T(t).] Integrating Eq. (|3}, re- 
placing min m [A m (T(i))] by A(T(t)), yields to 



T(t) 



pN 



k B \og(at + 1) 



o < t < r, 



(6) 



where a decreases exponentially with the system size N and 
is proportional to e, and T(T) is the temperature at which we 
want to study classical system. That is, if T is decreased as 
given by Eq. (|6jl, convergence to the desired state is guaran- 
teed. In the limit T(T) — > and logf > N > 1, we obtain 
T(t) ps (pN) /(kB log t) which agrees with the asymptotic 
convergence rate obtained in Ref. |9[ for SA. Such an agree- 
ment relies on the fact that the energy gap of H q (T) is also the 
energy gap of the transition matrix M q (T), which is known to 
determine the mixing time (or time required to reach thermal 
equilibrium) T M ps 0(1/A(T)). That is, in the SA scheme 
one never departs from equilibrium if the temperature is de- 
creased with the above convergence rate. Equivalently, in our 
context, the overlap between the adiabatically evolved quan- 
tum state and \tp(T(t))) is always close to 1. Note that Eq. (|6} 
holds even if the interactions in H are of long-range nature. 

QA has been proposed in Ref. Il3ll as an alternative method 
to reach the optimal (ground) state of a classical system with 
Ising-like interactions. Contrary to SA, the time-dependent 
quantum state in QA does not correspond, in general, to a 
thermal configuration of the original classical model. In this 
case, the quantum model Hamiltonian is given by H' (7) = 
H — 7 J2j a x> wnere 7 is decreased from a very large value, 
corresponding to T — ► 00, to 7 ps 0, corresponding to T ps 0. 
If 7 is slowly (adiabatically) changed, this method also al- 
lows us to reach the ground state of H. Similar techniques 
have been proposed to study the complexity of solving NP- 
complete problems, such as 3-SAT, using a quantum computer 
(QC) ITill . Numerical and analytical results show that, for cer- 
tain optimization problems, QA might enable a faster conver- 
gence rate to the optimal state than SA OH! El. Faster 
convergence of QA could be attributed to a decrease in the 
probability of driving the classical system to a local minima, 
as its dimension is effectively increased by one. Nevertheless, 
it has also been observed that in some cases 11711 QA performs 
similarly to SA. Note, however, that one could construct dif- 
ferent Hamiltonian paths to approach the optimal state. Each 
path yields to a particular convergence rate that has to be de- 
termined on a case by case basis. 

Using the classical-quantum mapping described above, the 
QA method can be extended to simulate classical statistical 



H q (l)=xYi 



lT,i <?l> having ^(7)) and \tp m (l)) 



as ground and excited states. Here, 7 is adiabatically de- 
creased from a very large value towards 7 ss x- in this way, 
the initial state J2[a] 1 1 "]) * s transformed into the desired state 
\t/j(T)). Notice that, from our viewpoint, QA differs from 
SA only by the choice of path used to reach the desired state. 
To successfully implement this annealing procedure, the rate 
at which 7 must be decreased is determined by the adiabatic 
condition, i.e. by the gap A(7) between the ground and first 
excited states of H q (^). This gap can be shown to satisfy 
A( 7 ) > 2-\/27riVe- iV '(l + c)- i V v = A( 7 ) HI, with 7 < c. 
Like the SA case, and for the worst-case scenario, the adia- 
batic condition H yields to 



7(t) ps [(2N - lXot)]- 1 ^*-!), < t < T, 



(7) 



where a depends on N, c, and e, and j(T) = x is determined 
by T. In the limit \ogt > N > 1, and 7(7") < 1, we ob- 
tain 7 (t) ps (2Nat)^ N \ If 1(^(7)1^(7)1^(7))! < 
xA m (7), with A m (7) the corresponding energy gap and x ps 
0(N q ), the coefficient 2N in Eq. |7} can then be replaced by 
N . In this manner, the con verg ence rate is in agreement with 
the result obtained in Ref. 11611 . Note, however, that this an- 
nealing schedule does not provide an advantage with respect 
to SA as 7 must be decreased to 7(T) = x> which is expo- 
nentially small in 1/T. 

The QA procedure to simulate T > can be directly im- 
plemented on a CC Hill . If the path-integral Monte Carlo 
method is chosen to simulate a d = 1 Ising-like model with 
nearest-neighbor interactions, H q (7) has to be mapped onto 
the 2-dimensional classical model, with energy functional 



E[a] 



N 



3ij(fi)a ik (T jk +^, t)J2Yl <7ika > 



k—l ij 



fe=l i=0 



(fc+1)- 

(8) 

Here, [a] = [cth, 021, • • • , <jnl] is one of the 2 N+L possi- 
ble spin configurations, and cr^ = ±1. The parameter L 
denotes the number of copies of the system in the extra di- 
mension (i.e., the Trotter discretization) and satisfies L ^> 1. 
The coupling constants Jij ((3) are defined via \ J2j e ^ Hj — 
ACS) + E y Jij(P)<r>i with Jij(fi 0) ps 0. The coeffi- 
cient (3 is given by the effective temperature of the quantum 
system and is not related with the temperature at which the 
classical system is studied. Therefore, $ 3> 1 and f3/L = St, 
with St being the time-slice of the discretization. The (fer- 
romagnetic) coupling between two adjacent copies is deter- 
mined by £(t) = \og[coth(Pxi(t) / L)]/2, and its magnitude 
increases as the transverse field 7(f) decreases to 7(T) = \, 
determined by T. In order to simulate more general classical 
systems at finite temperatures, the interactions appearing in 
Eq. (|8} must be modified accordingly. 

Note that the classical-quantum mapping can be extended 
and used to study any (finite-dimensional) classical system 
other than Ising-like models. In particular, it can be extended 
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to simulate s-spin classical systems (s > 1/2), where a vari- 
able <Tj = [— s, — s + 1, ■ ■ • , s — 1, s] is defined on each 
site j. In the case of QA, the ground state of Hg(j) = 
X J2j e ~ 1 Xj will determine the statistical proper- 
ties of the classical model when 7 — > x- The operators 
Xj £ bu(2s + l)j satisfy [Xj,S l z ] = 0, Vi 7^ j, and 
= -S|Xi, with X] = 1 'Here, e su(2)j is 
the angular momentum operator along the z-axis and deter- 
mines Hj. In matrix representation, X has l's in the anti- 
diagonal and O's otherwise. For example, in the s = 1 (three- 
state) Potts model O, E[a] = J £\ <K°j5 and ^ = 
J/2{E J [S£Si +1 (l+Si5|+ 1 )]-2[(5i) 2 + (5|+ 1 ) 2 ]}. There- 
fore, ffj = J/llSi^Sl + SiSi +1 ] and Xj- = 1 - (S^) 2 + 
[(S*^.) 2 + (S'i) 2 ]/2, defining the corresponding H° =1 ("f). 
The annealing schedule is again determined by adiabatically 
changing j(t) from a large value, where the initial state of the 
system is Yl[a] It 17 ])' t0 l^X) = X = e _/3p , where the final 
state of the system is e~^ H ' 2 53[cH I M)- 

It is important to stress that one can easily implement this 
extended QA (EQA) procedure by using current numerical 
methods. Since our analysis has only focused on the worst- 
case scenario, we would expect that for certain problems EQA 
should outperform SA JT3I1 . Moreover, one can always design 
other annealing procedures than the ones we have described. 
This can be done by constructing other quantum Hamiltoni- 
ans having ^(T)) as their ground state, and by introducing 
an extra interaction that is slowly changed to converge to the 
desired state. Depending on the path considered, it is expected 
a different behavior for the way that the relevant energy gap 
closes, and a different convergence rate as determined by the 
adiabatic condition. 

So far, we have considered that the lower bound in the 
gap is exponentially small in f3N, for SA, or TV log 7, for 
QA. One may wonder what the convergence rate for T(t) or 
7(f) is, when the gap can be bounded below by (/37V) -1 / 9 or 
(TV log 7) with q > independent of N and /3. In this 
case, integration of Eq. (0 yields to a convergence rate for SA 
of T(t) w 0[a/i 1 /te +1 )], with a a constant that depends on 
N and e. This is a much faster convergence rate than the one 
obtained in Eq. (|6). 

In this Letter, we have shown how to simulate the thermo- 
dynamic properties of an arbitrary classical model in d di- 
mensions by studying the ground state properties of a d di- 
mensional quantum system. This was achieved by an ex- 
act classical-quantum mapping. We have used the adiabatic 
theorem of quantum mechanics to analyze the convergence 
rate and resources required to reach the corresponding ground 
state. Our approach provides a unifying framework to address 
on an equal footing the well-known optimization methods of 
simulated and quantum annealing. These annealing proce- 
dures can be understood as two different evolution paths of 
the quantum system. It is remarkable, that the annealing rates 
obtained by using the adiabatic condition are in agreement 
with previous known results l9l Hill , which were obtained in 
the context of stochastic approaches such as path-integral or 



Green's function Monte Carlo. It is expected, however, that a 
QC will require less resources (e.g., quadratic speed-up) than 
a CC to solve these optimization problems. This issue will be 
addressed elsewhere. 
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